/* 
 * File:   main.c
 * Purpose:
 * Main Program
 * Author: abhi
 *
 * Created on September 21, 2011, 10:06 PM
 */

#include <stdio.h>
#include <string.h>
#include <stdlib.h>

#include "matrix.h"
#include "simplex.h"

/*
 * 
 */

/* Utility functions */
void getCostVector(int** costIntVector, int cols, int problemTypeFactor);
void getMatrixValues(int*** inputMatrix, int rows, int cols);
int frameMatrix(int*** inputMatrix, int rows, int cols, int** artificialVariables, int** basisVector, int* artVarNum);
rational getSolution(rational** solutionVector, rational*** rationalMatrix, rational** costVector, int** basisVector, int rows, int cols, int newCols, int artVarNum, int probTypeFactor);
void printFinalSolution(int rows, int cols, int newCols, int artVarNum, rational*** rationalMatrix, int** basisVector, rational** costVector, solution finalSolution, int probTypeFactor);
void printNote();

/* Main Program */
int main(int argc, char** argv) {
    
    int probTypeFactor;
    int rows,cols,newCols;
    int artVarNum;
    char probType[4];
    
    int **inputMatrix;
    int *artificialVariables;
    int *basisVector, *costIntVector;
    rational **ratMatrix, *costRatVector;
    solution finalSolution;
    
    /* get size of matrix */
    printf("Enter type of LP problem(max/min): ");
    scanf("%s",probType);
    printf("Enter number of variables: ");
    scanf("%d",&cols);
    printf("Enter number of constraints: ");
    scanf("%d",&rows);
    
    /* prepare cost Vector */
    if (strcmp(probType,"min")==0) {
        probTypeFactor = -1;
    } else {
        probTypeFactor = 1;
    }
    getCostVector(&costIntVector, cols, probTypeFactor); /* cols = variable columns */
    
    /* add column for value */
    cols++; /* cols = variable columns + value */
    
    /* get constraints as input */
    getMatrixValues(&inputMatrix, rows, cols);
    
    /* create matrix for simplex */
    newCols = frameMatrix(&inputMatrix, rows, cols, &artificialVariables, &basisVector,&artVarNum);
    
    printf("Initial");
    printIntMatrix(&inputMatrix, rows, newCols);
    
    /* extend cost vector */
    resizeIntVector(&costIntVector, (cols-1) , (newCols-1)); /* cols-1 and newCols-1 excludes value column */
    
    /* create rational cost vector */
    getRationalVector(&costRatVector, (newCols-1));
    convertIntToRationalVector(&costIntVector, &costRatVector, (newCols-1));
    
    /* create rational Matrix*/
    getRationalMatrix(&ratMatrix, rows, newCols);
    convertIntToRationalMatrix(&inputMatrix, &ratMatrix, rows, newCols);
    
    finalSolution = simplexProc( &ratMatrix, rows, newCols, &basisVector, artificialVariables, artVarNum, &costRatVector);
    
    printFinalSolution(rows, cols, newCols, artVarNum, &ratMatrix, &basisVector, &costRatVector, finalSolution, probTypeFactor);
    
    printNote();

    return (EXIT_SUCCESS);
}

void getCostVector(int** costIntVector, int cols, int problemTypeFactor) {
    
    int j;
            
    /* tmp Int vector */
    int* tmpIntVector;
    getIntVector(&tmpIntVector, cols);    
    
    printf("Enter cost function:\n");
    printf("Integral Coeff. of Variables from x1 to x%d\n",cols);
    
    for (j = 0; j < cols; j++) {
        printf("x%d :\t",(j+1));
        scanf("%d", &tmpIntVector[j]);
        tmpIntVector[j] = problemTypeFactor * tmpIntVector[j];
    }
    
    *costIntVector = tmpIntVector;
    
}

void getMatrixValues(int*** inputMatrix, int rows, int cols) {
    
    int i,j;
    char sign;
    
    /* create temporary matrix*/
    int** tmpInputMatrix;
    
    cols++; /* add for sign column */
    getIntMatrix(&tmpInputMatrix, rows, cols);
    
    /* Enter values */
    printf("Enter the coefficients of variables, Sign and Value for constraints below: \n");
    printf("(For Sign, please specify < and > for <= and >= respectively)\n");
    printf("Constraint\t");
    printf("Integral Coeff. of Variables from x1 to x%d, Sign, Value\n", (cols-2));

    for (i = 0; i < rows; i++) {
        printf("%d :", (i+1));
        for (j = 0; j < cols; j++) {
            
            if ( j == (cols-2) ) {
                printf("\t sign :\t");
                scanf(" %c", &sign);
                tmpInputMatrix[j][i] = sign;
            } else if ( j == (cols-1) ) {
                printf("\t value:\t");
                scanf("%d", &tmpInputMatrix[j][i]);
            } else {
                printf("\t x%d :\t",(j+1));
                scanf("%d", &tmpInputMatrix[j][i]);
            }            

        }
        printf("\n");
    }
    
    /* add artificial variables if required */
    *inputMatrix = tmpInputMatrix;
    
}

/*
 * Creates Matrix to be used for starting Simplex Method
 * 
 * returns Number of columns
 */
int frameMatrix(int*** inputMatrix, int rows, int cols, int** artificialVariables, int** basisVector, int* artVarNum) {
    
    int i,j,negative;
    int slackTotal=0,surplusTotal=0,equalTotal=0;
    int slackCount=0,surplusCount=0,equalCount=0,artificialCount=0,basisCount=0;
    int column=0;
    /* tmp Matrix */
    int **newMatrix;
    
    cols++; /* add for sign column */
    /* cols = variable columns + sign + value */
    
    /* determine if any column can be used in basis */    
    for (i = 0; i < rows; i++) {
        
        /* check if value column is negative */
        if ( (*inputMatrix)[cols-1][i] < 0 ) {
            negative = 1;
        } else {
            negative = 0;
        }
        
        for (j = 0; j < cols; j++) {
            
            if (negative) {
                if ( j == (cols-2)) {
                    switch ((*inputMatrix)[j][i]) {
                        case '<' : (*inputMatrix)[j][i] = '>'; break;
                        case '>' : (*inputMatrix)[j][i] = '<'; break;
                    }
                } else {
                    (*inputMatrix)[j][i] = (-1) * (*inputMatrix)[j][i];
                }
            }
            
            
            
        }
        
        switch ((*inputMatrix)[cols-2][i]) {
            case '<' : slackTotal++; break;
            case '=' : equalTotal++; break;
            case '>' : surplusTotal++; break;
        }
        
    }
    
    *artVarNum = equalTotal + surplusTotal;
    
    /* allocate to new Matrix */
    getIntMatrix(&newMatrix, rows, ((cols-1) + slackTotal + equalTotal + 2*surplusTotal ) );
    /* for slack and surplus variables and
     * artificial variables for equal and surplus variables
     */
    getIntVector(artificialVariables, (*artVarNum));
    getIntVector(basisVector, rows);
    
    
    /* slack and surplus variables given columns 1st 
     * then artificial variables assigned columns
     */
    
    /* (cols-2) accounts for last 2 columns of sign and value */
    for (i = 0; i < rows; i++) {
        
        /* check type of constraint */        
        switch ((*inputMatrix)[cols-2][i]) {
            case '<' : {
                /* column = cols-2+i; */
                column = cols-2+slackCount;
                newMatrix[column][i] = 1;
                slackCount++;
                (*basisVector)[basisCount] = column;
                basisCount++;
                break;
            }
            case '=' : {
                /* column = cols-2+slackTotal+surplusTotal+i; */
                column = cols-2+slackTotal+surplusTotal+equalCount;
                newMatrix[column][i] = 1;
                equalCount++;
                (*artificialVariables)[artificialCount] = column ;
                artificialCount++;
                (*basisVector)[basisCount] = column;
                basisCount++;
                break;
            }
            case '>' : {
                /* column = cols-2+rows+i; */
                column = cols-2+rows+surplusCount;
                newMatrix[cols-2+slackTotal+surplusCount][i] = -1;
                /* newMatrix[cols+slackTotal+surplusCount+equalTotal+i][i] = 1; */
                newMatrix[column][i] = 1;
                surplusCount++;
                (*artificialVariables)[artificialCount] = column ;
                artificialCount++;
                (*basisVector)[basisCount] = column;
                basisCount++;
                break;
            }
        }
        
    }
    
    cols = cols - 2; /* remove sign and value columns from cols */
    /* cols = variable columns */
    
    /* copy matrix to newMatrix */
    /* copy the variable columns */
    copyIntMatrix(inputMatrix, &newMatrix, 0, 0, 0, 0, rows, cols);
    /* copy sign and value columns */
    copyIntMatrix(inputMatrix, &newMatrix, 0, (cols+1), 0, (cols + slackTotal + equalTotal + 2*surplusTotal), rows, 1);
    
    /* replace input Matrix with new Matrix */
    (*inputMatrix) = newMatrix;
    
    return (cols + slackTotal + equalTotal + 2*surplusTotal + 1); /* 1 for value column */
    
}

rational getSolution(rational** solutionVector, rational*** rationalMatrix, rational** costVector, int** basisVector, int rows, int cols, int newCols, int artVarNum, int probTypeFactor) {
    
    int i;
    int sum;
    rational product;
    rational* newCostVector;
    
    /* truncate cost vector and restore according to LP Problem type */
    getRationalVector(&newCostVector, (newCols-1));
    copyRationalVector(costVector, &newCostVector, (newCols-1));
    resizeRationalVector(&newCostVector, (newCols-1), (cols-1));
    resizeRationalVector(&newCostVector, (cols-1), (newCols-1));
    getRationalScalarRationalVectorPdt(itor(probTypeFactor), &newCostVector, (newCols-1), &newCostVector);
    
    /*printRationalVector(&newCostVector, (newCols-1));
    printRationalMatrix(rationalMatrix, rows, (newCols-artVarNum));*/
    
    for(i=0; i<rows; i++) {
        (*solutionVector)[(*basisVector)[i]].numerator = (*rationalMatrix)[(newCols-artVarNum-1)][i].numerator;
        (*solutionVector)[(*basisVector)[i]].denominator = (*rationalMatrix)[(newCols-artVarNum-1)][i].denominator;
    } /* length of basisVector = # rows in rational Matrix */
    
    /*printRationalVector(solutionVector, (newCols-1));*/
    
    /* compute cost */
    product = itor(0);
    for(i=0; i<(newCols-1); i++) {
        product = rAdd(product, rMul(newCostVector[i], (*solutionVector)[i]));
    }
    
    return product;
    
}

/*
 * Prints Final Solution
 */
void printFinalSolution(int rows, int cols, int newCols, int artVarNum, rational*** rationalMatrix, int** basisVector, rational** costVector, solution finalSolution, int probTypeFactor) {
    
    int i;
    rational solutionCost;
    rational* solutionVector;
    
    for(i=0;i<15;i++) { printf("-"); }
    printf("\n\n");
    printf("Final Solution:\n\n");
    
    switch (finalSolution) {
        case OPTIMAL: printf("Solution is Optimal"); break;
        case UNBOUNDED: printf("Solution is Unbounded"); break;
        case INFEASIBLE: printf("Problem is Infeasible");
    }        
    printf("\n\n");
    
    if (finalSolution == OPTIMAL) {
        
        /* create solution vector */
        getRationalVector(&solutionVector, (newCols-1)); /* value column excluded */
        
        solutionCost = getSolution(&solutionVector, rationalMatrix, costVector, basisVector, rows, cols, newCols, artVarNum, probTypeFactor);
        /* solutionCost = basisCostVectorPdt(rows, basisVector, costVector, rationalMatrix, (newCols-1)); */ /* valueColumn index = (newCols-1) */ 
    
        printf("\nSolution Cost = %d/%d \n",solutionCost.numerator,solutionCost.denominator);
        printf("Solution is :\n");

        /* for all variables given */
        for(i=0; i<(cols-1); i++) {
            printf("x%d = %d/%d\n",(i+1),solutionVector[i].numerator,solutionVector[i].denominator);
        }
        
    }
    
    printf("\n");
    
}

/*
 * Prints Explanatory Note
 */
void printNote() {
    
    printf("Note: In the above matrices:\n");
    printf("1. Rational Numbers are represented in the form (p/q) where q!= 0 \n");
    printf("2. Zero is represented as 0/1 \n");
    
}